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The Hopfield network is reconstructed as an inverse Ising problem by passing messages. The 
applied susceptibility propagation algorithm is shown to improve significantly on other mean-field- 
type methods and extends well into the low temperature region. However, this iterative algorithm is 
limited by the nature of the supplied data. Its performance deteriorates as the data becomes highly 
magnetized, and this method finally fails in the presence of the frozen type data where at least two 
of its magnetizations are equal to one in absolute value. On the other hand, a threshold behavior is 
observed for the susceptibility propagation algorithm and the transition from good reconstruction 
o : to poor one becomes sharper as the network size increases. 
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I. INTRODUCTION 

Message passing algorithms have important applications in various contexts ranging from random constraint satis- 
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faction problems [l|, supervised learning [2| to information theory [H, [j]. In recent years, the idea of message passing 
was introduced into the study of the inverse Ising problem [5|-l7|. Given observed magnetizations (mean activities) 
and pairwise correlations {m^CV,}, the inverse Ising problem aims at finding the underlying parameters {hi, Jij}, 
local fields and coupling constants to describe the statistics of the experimental data which can be collected either 
in real experiments (e.g., microarray measurements in gene expression experiments or multi-electrode recordings in a 
I • neuronal population) or in Monte Carlo simulations. The active research of inverse Ising problem is mainly motivated 
by the observation of correlated activity in the retinal network [8l-[l0| , the cortical network [l2j and other biologi- 
cal networks [H3rfi~5l |. Encouragingly, the pairwise Ising model, as a least structured model, was shown to be capable 
of capturing most of the correlation structure of the network activity [ill, EH, EH E3 ■ Based on the Ising model, 
computationally efficient inverse algorithms were proposed to analyze multi-electrode recordings in the salamander 
retina (l8j and to identify correlation between amino acid positions in interacting proteins (l9| . The pairwise Ising 
model only requires 0(N 2 ) parameters to describe the original distribution and is thus attractive for dimensional 
reduction in modeling vast amounts of biological data (for a general statistical physics analysis of the pairwise Ising 
' model on the inverse Ising problem , see, e.g., Ref. 

For large system, the inverse Ising problem is known to be a hard computational problem [2iT[23j . Various approx- 
imate schemes were proposed to tackle this problem. As one of these approximations, message passing strategy looks 
, promising and the susceptibility propagation (SusProp) algorithm has been derived to infer couplings of Sherrington- 
Kirkpatrick (SK) model @. SusProp solves a closed set of equations by iteration, passing messages along the directed 
edges of the factor graph representation [24[ of the problem. The new message is computed only based on the in- 
■ coming messages. This locally updating feature makes SusProp fully distributed and amenable to parallelization. 
Correlation information of any two variables provided by SusProp can be used not only to decimation procedure 
in solving some hard constraint satisfaction problems (25[ but also to the network reconstruction [(| 0, Ell- Depen- 
dence of the performance of SusProp on the quality of the originally observed data has been studied in Ref. 
At high temperature, the quality of reconstruction is constrained by the implementation precision of the algorithm 
and the random noise embedded in the supplied data. Statistical errors presented in the Monte Carlo noisy data 
could also have detrimental effects on the reconstruction performance. Aurell et al. in a recent work [6] studied the 
dynamical behavior of SusProp on inferring couplings from synthetic data of SK model. They found that, at the 
low temperature (T < 4.0), the algorithm doesn't converge typically with diverging inferred couplings, however, by 
introducing a stopping criterion, the threshold could be pushed to lower value. A transition from reconstructible to 
non-reconstructible phase for SusProp was also observed in their numerical simulations. High absolute magnetization 
was claimed to have negative effects on the performance of the algorithm. All aforementioned investigations were 
restricted to the SK model. Mean field schemes based on inversion of correlation matrix have been recently tested on 
Hopfield networks [26]. Regarding these mean field schemes, the simple ones are naive mean field (nMF) method and 
independent pair (ind) approximation, and the more advanced ones inversion of Thouless-Anderson-Palmer (TAP) 
as well as Sessak-Monasson (SM) approximation (for details, see Refs. [26|, |27| . a brief description is also given in 
Appendix [A| . It was shown that all mean field schemes fail to extract interactions within a desired accuracy in the 
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retrieval phase. Zhang et al. (28j applied recently belief propagation plus an auxiliary updating external field to infer 
couplings of the sparse Hopfield network when the system settles in the retrieval phase. They showed that inference 
error with sampling from single basin of the stored pattern is much similar while error with sampling from multiple 
basins is drastically reduced. 

In the present work we will examine the reconstruction performance of SusProp on both the fully connected and 
sparse Hopfield networks and discuss the limitation of this message passing algorithm. Improvements over other 
existing mean field methods are reported and a threshold behavior relative to SusProp is observed in our simulations 
on single instances. The rest of this paper is organized as follows. The definition of the Hopfield network is given in 
Sec. |lT] followed by the detailed demonstration of SusProp in Sec. Mil In Sec. IIV1 we report improvements achieved 
by SusProp and discuss its limitation and threshold behavior. Conclusions and future perspectives are devoted to 
Sec. El 



II. HOPFIELD NETWORKS 

The Hopfield model was proposed to mimic the memory and recall functions of real neuronal networks p9l I30I ] . 
It yields rich statistical physics properties and is moreover a simple model to assess the efficiency of various inverse 
algorithms [26| . which may have some implications for neuroscience. The equilibrium properties of the Hopfield 
network are governed by the following Hamiltonian: 

H = -^2,JijVi<Tj (1) 

i<j 

where Oi describes the state of each neuron in the network. <n = +1 indicates the spiking of neuron i and <jj = — 1 
the silence of the neuron. Coupling Jy is constructed according to the Hebb's rule: 

(2) 

where {£" } taking ±1 with equal probability are P stored random patterns. The ratio of the number of stored patterns 
to the network size N is termed the memory load of the network, i.e., P — aN. In the fully connected network, each 
neuron is connected to all the other neurons and no self-interactions are assumed. The mean field behavior of the 
fully connected Hopfield model has been thoroughly studied in Refs. [3l|, |32|]. When a < 0.138, paramagnetic, spin 
glass and metastable ferromagnetic retrieval phases appear in order as the temperature decreases. The retrieval phase 
becomes stable at low temperatures if a < 0.051. Replica symmetry breaking occurs for the retrieval phase only at 
very low temperatures. However, the replica symmetry solution for the spin glass phase are unstable in the entire 
spin glass phase [33j | . For the sparse network, the coupling or interaction is constructed as 

■', y£c^ (3) 

where I is the mean degree of each neuron. In the thermodynamic limit, P scales as P = al where a is the memory 
load. No self-interactions are also assumed and the connectivity Uj follows the distribution: 

P(hj) = (i - ^) Kki) + jf—^i - 1) (4) 

Mean field properties of the sparse Hopfield network have been discussed within replica symmetric approximation in 
Refs. [3~il l35l|. Three phases (paramagnetic, retrieval and spin glass phases) have also been observed in this sparsely 
connected Hopfield network with arbitrary finite I. For large I (e.g., I = 10), the phase diagram resembles closely that 
of extremely diluted case [3H where the transition line between paramagnetic and retrieval phase is T = 1 for 
a < 1 and that between paramagnetic and spin glass phase T = \fa for a > 1. The spin glass/retrieval transition 
occurs at a = 1. 

We simulate the Hopfield network using Glauber dynamics plus simulated annealing techniques to collect enough 
experimental data (totally 10 samplings): = (ci) data , Cy = ("iCi^ata — mirrij where (• • ■ ) data denotes the aver- 
age over the collected data (simulation details are given in Appendix |B|) . then we estimate the couplings { Jy} between 
neurons based on these measured magnetizations and two-point connected correlations, such that the resulting Ising 



distribution P\ s i n ^ cx exp 



■ J*-OiOj is able to provide an accurate description of the statistics of the experimental 
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data. Note that J*j is the inferred value and has been scaled by the inverse temperature j3. The external field hi 
is zero for all neurons in the current Hopficld model. In other cases [HI, [III [IH, the external field can be used to 
represent the preferred direction of cr, . 



III. INFERRING COUPLINGS BY PASSING MESSAGES 



Neurons in the Hopficld network usually interact with each other to yield collective behavior at the network level. 
{Cij} measure the tendency for each pair of neurons to spike cooperatively. Together with the information about mean 
firing rates {m,}, they can be used as inputs to SusProp for the network reconstruction. SusProp passes messages 
along the directed edges of the network by iterative updating. To iterate SusProp, two kinds of messages are needed. 
We first define the cavity magnetization rrii^j as the message propagating from neuron i to neuron j, then define the 

other kind of message, namely the cavity susceptibility gi->j,k = q^ 3 where hi->j is termed cavity field of neuron i 
in the absence of neuron j and hk is the local perturbation. The SusProp then reads as follows: 



m, 



Ji— >j,k 



rrii — mj->i tanh Jy 



1 — miirij^i tanh J,j 
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(5a) 
(5b) 

(5c) 
(5d) 



where di\j denotes neighbors of neuron i except j, Sik is the Kronecker delta function and e serves as a damping 
factor. The damping factor allows the new to memorize a given fraction (1 — e) of Jy computed at the last step 
(denoted as J°j d ), which helps SusProp converge to a fixed point if the temperature is not very low although the 
convergence is slowed down. In our current simulations, we employ very small e of order from 0(1CP 2 ) to 0(1O~ 4 ). 
Detailed derivation of SusProp is given in Appendix [C] For other discussions of this algorithm, we refer the reader to 
previous works 0-0, HI- The SusProp algorithm is able to estimate the correlation between any two variables even 
if they are not directly linked in the network [25j] and this information could be used further to update couplings. 
Furthermore, the memory term e ensures the update of Jij towards its true value step by step when the temperature 
is not very low. For the fully connected network, SusProp has the complexity of 0(N 3 ). 

To assess the reconstruction performance of SusProp, we define intuitively the inference error as 



A 



W* rtruc\2 
ii J a ) 



N(N-l) ^ 



1/2 



(6) 



where J*j is the inferred value of the coupling and Jfj Ue the original coupling constructed according to the Hebb's rule 
Eq. ([2]) or Eq. ([3]). To run SusProp, we initially set all Jij to be zero, and randomly initialize for every edge of the 
network the message rrii^j € [—1.0, 1.0] and g^j.k — if i ^ k and 1.0 otherwise. Then SusProp is iterated according 
to Eq. ([5]) until the inferred couplings converge or the preset maximal number of iterations %nax is saturated. In 
our simulations, we adopt convergence criterion 77 = 10~ 4 , i.e., convergence of SusProp is identified once all updated 
couplings have converged within precision rj. 



IV. RECONSTRUCTION PERFORMANCES AND DISCUSSIONS 



To avoid the expensive computational cost, we assess the reconstruction performance of SusProp only on small 
size networks. The phase diagram of the Hopfield model has been derived for the fully connected case [HI, [32| and 
the finite connectivity case [34| . In the finite size system, we distinguish different phases by two order parameters; 
one is the overlap between the network configuration and the /zth pattern m M = -h X)i=i £f m » an( i the other is the 



4 




FIG. 1: (Color online) Comparison of reconstruction performances of SusProp and other existing mean field schemes in the 
fully connected Hopfield network. The inference error is plotted against the memory load with T = 0.6, N = 100. Lines are 
guides to the eye. Each point is an average over five random samples and error bars are also shown. 




FIG. 2: (Color online) The inference error versus temperature for the fully connected network with N = 100. Lines are guides 
to the eye. Each point represents an average over five random samples and error bars are also shown, (a) Comparison of 
reconstruction performances of SusProp and other mean field schemes for a — 0.03. (b) Inference error versus temperature for 
SusProp with different memory loads. 



mean-squared magnetization q — X)i=i m i- To implement SusProp, 7~ m ax is set to be 2000 and we need to adopt 
an appropriate damping factor to prevent the absolute updated tanh Jij from being larger than one. 

Comparison of reconstruction performances of SusProp and other mean field schemes is shown in Fig. [TJ SusProp 
turns out to be the most efficient, reducing the inference error by a significant amount. Moreover, it seems to be 
less sensitive to the memory load compared with other mean field methods. Fig. [2] reports the inference error as 
a function of temperature for various reconstruction algorithms. SusProp operates fairly accurately and extends 
the reconstructible region well into a much lower temperature down to 0.4. However, when the system settles in 
low temperatures (T < 0.6), SusProp first suffers highly magnetized data and its performance gets worse with non- 
convergence, which may be remedied by adopting much smaller e and larger %nax or by introducing a stopping 
criterion jg]. When the temperature approaches lower values (e.g., below 0.4 in Fig. [5] (a)), the collected data will 
then become frozen, i.e., at least two of its magnetizations equal one in absolute value, as a result, SusProp yields 
diverging couplings and fails to reconstruct the network. Actually, in the high absolute magnetization case, both |CV, | 
and Imi-^l (and |mj_H|) are very close to one on some edges < ij >. We then rewrite CV, as CV, = sgn(CV,)(I — £y ), 
and similarly mi^jrrij^i = sgn(rrii-^jmj^i)(l — where sgn(-) is a sign function; Sij and Eij^i are small positive 
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FIG. 3: (Color online) Inference performances for the sparse Hopfield network with a = 0.6, 1 = 5, N = 100. Lines are guides 
to the eye. Each point is an average over five random samples and error bars are also shown. Efficiencies of SusProp and other 
mean field schemes are compared, so are the reconstruction performances of SusProp with and without prior knowledge on the 
sparseness of the network. 




FIG. 4: (Color online) Fraction of good reconstruction of fully connected Hopfield networks as a function of inverse temperature 
/3 for SusProp. Different network sizes are considered with the same memory load a = 0.1. Twenty random samples are 
simulated and a good reconstruction is identified when A < 0.030. 



values compared to one. Then one can readily recast tanh Jy according to Eq. (JS|) as 

tanh J, = S ""' r '-' S I _ 

for Cijmi^jmj^i > and | tanh Jy| = 1.0 for Cijmi^jmj^i < 0. If Cijmi^jmj^i > and neither of Ey- and £ijji 
vanishes, the estimated remains finite and SusProp does work. In other cases (e.g., |Cy| = 1.0 and \mi_+jmj^i\ ^ 
1.0), Jij suffers divergence and SusProp is unable to infer the couplings. In this situation, at least two of the supplied 
magnetizations are equal to one in absolute value, e.g., \m,i\ = 1.0 and \rrij\ = 1.0, under the update rule Eq. ([S]), 
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nii^j — sgn(mi) and rrij^i = sgn(m 3 -). These highly polarized messages will then spread out in the network, which 
yields an infinite value for as explained above. A simple physical interpretation is that, in this case, CV,- ~ rriivn^ 
and this implies that some neurons in the network tend to behave independently of other neurons and thus SusProp 
couldn't get all information about correlations of the network, which leads to the failure of reconstruction. In the 
current Hopfield model, the frozen type data does appear for small enough a and T where the system gets trapped by 
one of stable or metastable memory states (26l |. Increa sing a but maintaining low T, many spurious minima show up 
and Monte Carlo sampling becomes very difficult [IH, l38| which produces fairly noisy collected data. On the other 
hand, the current SusProp hasn't taken into account the complex structure of the phase space, therefore it remains a 
non-trivial issue for SusProp to deal with this more involved case. The reconstruction error against temperature is also 
shown with respect to different memory loads for SusProp in Fig.[5](b). In the high temperature region (T > 0.6), the 
smaller the memory load is, the more precisely SusProp reconstructs the Hopfield network. As temperature decreases 
further, SusProp becomes less precise for all memory loads while still maintaining a relatively small error. 

SusProp is applied to reconstruct the fully connected Hopfield network, whereas, it shows a surprisingly good 
performance. If the network is sparse and locally treelike, SusProp is believed to be fast and able to give a precise 
estimation. To test its efficiency on reconstructing the sparse network, we compare performances of SusProp with 
those of other mean field schemes in Fig. [3l It is clearly shown that SusProp performs exceedingly well particularly 
in the low temperature region (down to T — 0.5) and exhibits less sample-to-sample fluctuations. If we have a prior 
knowledge of the sparseness of the network, i.e., we know a priori the connectivity pattern of the network, the inference 
error could be reduced substantially. In this case, we only infer the strength of interaction between neurons which 
are really connected. In fact, when the system is presented at the high temperature, one can reconstruct the network 
using an appropriate cutoff since the estimated couplings between unconnected neurons are very small compared to 
those between really connected neurons. 

In Fig. 21 we report the fraction of good reconstruction versus temperature for different network sizes at fixed 
memory load a — 0.1. A good reconstruction is identified when A < 0.030 and SusProp converges within T ma x- As 
increasing temperature, a threshold behavior is observed and the transition becomes sharper with growing network 
size. Given large enough N, the probability that SusProp gives good reconstruction tends to be zero when the 
temperature is below the critical value, and at a high enough temperature, SusProp succeeds in reconstructing the 
network in all instances. The critical temperature is estimated to be about 0.6. For a more precise estimation, more 
samples and larger network size are required. The threshold behavior of SusProp on network reconstruction has also 
been observed in SK model [6j. 



V. CONCLUSION AND OUTLOOK 



SusProp solves the inverse Ising problem by iteratively updating messages along the directed edges of the network, 
and is shown in the present work to outperform all other mcan-ficld-type schemes and extend the reconstructible region 
into lower temperatures. We also study the fraction of good reconstruction as a function of temperature and a threshold 
behavior is observed. The transition from good reconstruction to poor one becomes sharp as increasing network size 
and the critical temperature is estimated to be about 0.6. SusProp is also amazingly efficient for reconstructing sparse 
Hopfield network and its performance can be further improved by introducing a prior knowledge about the sparseness 
of the network. The sparse case is more relevant in modeling real biological data than its dense counterpart. We hope 
our analysis of SusProp on the Hopfield network reconstruction can be extended to the more biologically relevant 
cases. 

At high temperatures, the performance of SusProp is believed to be limited by the quality of the supplied data Q. 
Once the system is presented at the low temperature, the efficiency of SusProp is also determined by the nature of the 
input data. In the presence of highly magnetized data, the reconstruction performance of SusProp gets deteriorated 
with a high inference error. Furthermore, the frozen type data makes updated couplings diverge and any value of 
damping factor can not overcome this hurdle, reminiscent of the fact that the frozen phase in random constraint 
satisfaction problems is most difficult for any known algorithm [38], |39| . 

For large enough a but low enough T, the system enters the spin glass phase where Glauber dynamics is easily 
trapped by one of the spurious minima correlated or uncorrelated with the stored patterns. SusProp fails to extract 
couplings precisely in this region since it does not take into account the complex structure of the phase space. In fact, at 
finite temperatures, the support of cavity field distributions becomes real- valued and a sampling procedure is required. 
On the other hand, other values of Parisi parameter (smaller than the optimal value associated with the ground states) 
also carry physical information and can be used to describe the metastable states which trap Glauber dynamics [40j . 
However, searching for an optimal Parisi parameter (also known as replica symmetry breaking parameter) is also a 
time consuming task for the network reconstruction [4l| . Rather, provided that Glauber dynamics gets stuck in some 
metastable state for a very long time, how much information can be extracted from this state about the couplings of 



7 



the network remains an important issue for future study. 
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Appendix A: Mean field schemes for network reconstruction 
1. Naive Mean-Field Method 



The naive mean field theory gives rrii = tanh [hi + X^fe^i JikmkJ where hi is the external field and 
magnetization. Using the fluctuation-response relation, 



dhj 

one obtains the nMF prediction of couplings, 



Ca = = (1 - mj 



k^i 



where Py = (1 — rn 2 )<5y . 



j^ = (p-%-(c-i) i 



2. Independent-Pair Approximation 



(o-j) the 



(Al) 



(A2) 



In this approximate scheme, each pair of neurons are independent of other neurons of the system, i.e., their joint 



probability P(<7i, <jj) oc exp 



h\ <Ji + hj <Jj + JijtJiGj 



where h^(hj ) is the local field neuron i(J) feels when neuron 



j(i) is removed from the system. Then the ind prediction is given by 



ind 



log 



(l-C;.) 2 -K-m,) 2 



(A3) 



where Cj = dj + rrumj. 



3. Sessak-Monasson Approximation 



The SM prediction of couplings is derived based on a perturbative expansion in the correlations [42J and it can be 
formulated as 



J 



SM 



J, 



nMF 



J 



ind 



(l-m?)(l-m?)-C? 



(A4) 



4. Inversion of TAP Equations 

The usual TAP equation reads hi = tanh -1 to, — Y^j^i Jij m j + m i Yj^i ^ijO- ~ m j) 0- Differentiating the field 
hi with respect to the magnetization rrij, one readily obtains the TAP prediction equation, 



dhj 
drrii 



(A5) 
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Appendix B: Simulation details 



The rule for Glauber dynamics can be generally expressed as P(<Ji — > — 0i) = 1+exp ^ A ^ ^ where ATii is the energy 
change due to such a flip. For the current Hopfield model, the dynamics rule is recast into 

P{(Ti -> -a,) = - [1 - a, tanh Ai] (Bl) 

where /3 is the inverse temperature and hi = ■^ij (J j i s the local field acting on crj. 

In our numerical simulations, we update the state of each neuron according to Eq. IBll in a randomly asynchronous 
manner. We define a Glauber dynamics step as N proposed flips. Introducing simulated annealing strategy, we set 
the initial temperature to be 1.0 and the cooling rate 0.005. At each intermediate temperature, we run 10 4 Glauber 
dynamics steps. When the temperature is decreased to the desired one, we run another 2 x 10 6 steps to calculate 
magnetizations and correlations. We sample the state of the network every 200 steps. For high temperatures (> 1.0), 
we run totally 4 x 10 6 steps, among which the first 2 x 10 6 steps are run for the system to reach the equilibrium state 
and the other 2 x 10 6 steps for calculating magnetizations and correlations. 



Appendix C: Derivation of SusProp update rules 

To derive Eq. [SJ we first derive the susceptibility propagation equations for general if-body interaction problem 
(K = 2 for the Hopfield model). Using the factor graph representation [24], we denote a as the function node 
representing the constraint imposed on a subset of spins (jQ a (da denotes neighbors of the function node a) , and i as 
variable node representing the spin on the factor graph. The belief propagation is then formulated as 0] 



mi 



>a = tanh/ii_j. a = tanh Ub-rt — tanh [tanh '(m,) — w a -n] 

i b£di\a 



tanh JaUjeda 



\i m j^a 



1 - to* tanh JaUjeda 
tanh Ufc_>.i = tanh Jf, j | mj^b 



\i tnj^a 



(Cla) 



(Clb) 



jedb\i 



where h^ a is the cavity field (correspondingly m^a i s the cavity magnetization) acting on spin <7j in the absence of 
a; the cavity bias when i is involved in b only, and hi^ a , as well as Jj, have been rescaled by (3. 

We define cavity susceptibility g^ a .k = . From the belief propagation equations, one readily gets the update 

rule for g^ a ,k- 



9i^a,k — $ik + ^ 



du b - 



b£di\a 



dh k 



= S ik + 



tanh Jfc 



b€di\a 1 - (tanh Jb Y[jedb\i m 3^b) j£db\ t 



E 



Y[ m n-ib 
nEdb\i,j 



g 3 ^ b ,k{l - rn]^ b ) 



(C2) 



Using the identity [25j: 



du a -yi ({hj^ a }j£da\i) 
dh n ^a 



(a n (jj) - (a„ ) {<Ji 



1 - (<*} 



v 2 



(C3) 



where i,n € da and (•••) denotes the average under the joint probability distribution Prob(crg a ) which can be 
computed from the belief propagation equations, one can re-express the correlations Cij through the fluctuation- 
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response relation: 



_drrii _ t , 



dhj 



(1 - mi) 



E 



bedi 3 



(1 - 



(1 - m\) 



du a —>i dh n ^ a 
dh n ^ a dhj 



(C4) 



For the simple case, the Hopfield network, J Q = J^, and the constraint a involves only two neurons, say i 

and j. To obtain the update rule for Jy, we compute Cij = {&i<Jj) directly by assuming Prob(<7i , a 3 ) oc 
exp (JyiTjiTj + hi^jai + h^iCTj), finally we get 



tanh Ji. 



Tfli^yjTTlj— 

1 - ( ', ,///, . ,;// , 



(C5) 



where Cj,- can be evaluated from Eqs. IC3l and lC41 i.e., 



From Eqs. ICll and IC21 the update rules for and gi^^k are finally obtained as follows, 



(C6) 



rrii — rrij^i tanh Jy 



1 — miVdj^i tanh J^. 



1 — raz 



iedi\j 



(m n ^i tanh J„ 



■ tanh JniSn—H,k 



(C7a) 
(C7b) 



Introducing additionally a damping factor e, Eqs. IC51 IC6I and IC7I are the very SusProp we have presented in Sec. IIIII 
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